Fish_Stuff
Callie Stephenson
2024-03-28
knitr::opts_chunk$set(echo = FALSE, message = FALSE, warning = FALSE)TIDY DATA:
First chunk of code subsets to only using 8 blocks per survey location, 4 from an ‘on’ and 4 from an ‘adjacent’. I also remove pins where camera malfunctions made this not possible (Pin 13 and Pin 2).
Removed unknowns, cryptic species, etc:
Since we could count some species as distinct individuals by life stages, I have combined those counts together here:
And built a little code to double check that it went correctly (should create two empty data frames)
I then deal with extreme schooling events following methods from the Donovan regimes paper: “Additional methodology was developed for dealing with outliers in the fish data, accounting for extreme observations of schooling species. Extreme observations in the database were defined by calculating the upper 99.9% of all individual observations (e.g. one species, size and count on an individual transect), resulting in 26 observations out of over 0.5 million, comprised of 11 species. The distribution of individual counts in the entire database for those 11 species was then used to identify observations that fell above the 99.0% quantile of counts for each species individually. These observations were adjusted to the 99.0% quantile for analysis.”
I want to make sure blocks are summarized appropriately, so I renamed them:
and ending by renaming the dataframe as just fish_tidy for future use: this allows me to make extra tidy data adjustments without breaking the dataframe name later on in the code!
Resource Fish Only
I’m pulling a list of fished genus from the paper “Perceptions and responses of Pacific Island fishers to changing coral reefs” by Rassweiler et al., 2021. I am editing fish_tidy to split genus and species, creating a column based off whether the genus is on the list for Rassweiler et al, and then making that into a binary variable ‘Resource’.
Also in this paper, it shows that “Naso, Chlorurus and Scarus collectively composed the bulk of the fishable biomass on the reef (48–66%) and a roughly similar total proportion of the catch (43–65%).” I have made a variable ‘most_targeted’ for these genus.
Wide data for ordination:
I’ve also made a wide version of the data frame for ordination:
Then used that one to re-make the long data frame with a presence column instead of abundance
Data Vis
Linda from Johanssen Lab reccommended to make a plot that is ordered by the most abundant fish near the seep and look for any patterns, so I went and made this plot a number of ways. First, I kept all the species and looked a just the average abundance at all pins, to see if the communities looked relatively similar in this way.
All species average abundance across gradient
This is too busy to read, and no major patterns jump out. I subset this to just the top 20 species along the gradient, and again the order of the pins does not have biological relevance.
Most common 20 species average abundance across gradient
So, now I will order the pins by those most common at pin 14
This isn’t showing the most obvious patterns, but there were only a few species of fish that swam in the sand near the seep. So, I’ll reorder the pins by the distance to the seep
Most common species average abundance at each pin by distance to seep:
To see if maybe a nutrient parameter would be good descriptor, I do this also by the low tide mean silicate
Most common species average abundance across the gradient of low tide mean silicate
Most common species average abundance across the gradient of nutrient delivery
Nutrient delivery quantified as a pc axis I got from the low tide mean of major nutrients across the gradient.This axis captures 72.87298% of the variation in mean low tide TA, Phosphate, Silicate, and N + N. Higer numbers on this axis indicate higher values of these variables.
Alpha Diversity:
Calculating richness against gradient metrics:
Calculating diversity against gradient metrics:
NMDS Non-Dimensional Ordination
Code from Linda:
I’ve been going back and forth on the transformation here, and am currently going with hellinger instead of vegan’s default, the wisconsin. This is because, from my limited understanding, the hellinger better handles cases where species are not present in multiple sites (known as double zeros).
## Run 0 stress 0.2193151
## Run 1 stress 0.2286818
## Run 2 stress 0.2303236
## Run 3 stress 0.2154216
## ... New best solution
## ... Procrustes: rmse 0.02448041 max resid 0.1813317
## Run 4 stress 0.2182822
## Run 5 stress 0.2243652
## Run 6 stress 0.2279111
## Run 7 stress 0.2237805
## Run 8 stress 0.2151307
## ... New best solution
## ... Procrustes: rmse 0.02295824 max resid 0.1673727
## Run 9 stress 0.2259491
## Run 10 stress 0.2225274
## Run 11 stress 0.2273464
## Run 12 stress 0.2195989
## Run 13 stress 0.2236718
## Run 14 stress 0.2176085
## Run 15 stress 0.2259154
## Run 16 stress 0.2303317
## Run 17 stress 0.2248568
## Run 18 stress 0.2193465
## Run 19 stress 0.2207669
## Run 20 stress 0.2211457
## Run 21 stress 0.2236639
## Run 22 stress 0.2240872
## Run 23 stress 0.2161173
## Run 24 stress 0.2212809
## Run 25 stress 0.2155593
## ... Procrustes: rmse 0.01960847 max resid 0.1704885
## Run 26 stress 0.2247331
## Run 27 stress 0.2228412
## Run 28 stress 0.2290937
## Run 29 stress 0.2251805
## Run 30 stress 0.2278835
## Run 31 stress 0.2319124
## Run 32 stress 0.2320552
## Run 33 stress 0.221125
## Run 34 stress 0.22673
## Run 35 stress 0.2267868
## Run 36 stress 0.2232394
## Run 37 stress 0.2161406
## Run 38 stress 0.2270662
## Run 39 stress 0.2209657
## Run 40 stress 0.2237046
## Run 41 stress 0.2197634
## Run 42 stress 0.2229924
## Run 43 stress 0.22823
## Run 44 stress 0.2272547
## Run 45 stress 0.2162016
## Run 46 stress 0.2280137
## Run 47 stress 0.231659
## Run 48 stress 0.2160515
## Run 49 stress 0.2237892
## Run 50 stress 0.2231962
## Run 51 stress 0.2152797
## ... Procrustes: rmse 0.009919574 max resid 0.08774457
## Run 52 stress 0.2255716
## Run 53 stress 0.2160574
## Run 54 stress 0.2190227
## Run 55 stress 0.2193456
## Run 56 stress 0.2296862
## Run 57 stress 0.2254759
## Run 58 stress 0.2251675
## Run 59 stress 0.2225763
## Run 60 stress 0.2155459
## ... Procrustes: rmse 0.02072545 max resid 0.1699998
## Run 61 stress 0.2284182
## Run 62 stress 0.2201268
## Run 63 stress 0.2230935
## Run 64 stress 0.2182881
## Run 65 stress 0.2220194
## Run 66 stress 0.2311721
## Run 67 stress 0.221164
## Run 68 stress 0.2187999
## Run 69 stress 0.2290145
## Run 70 stress 0.2206629
## Run 71 stress 0.2349865
## Run 72 stress 0.2255399
## Run 73 stress 0.2258342
## Run 74 stress 0.2253595
## Run 75 stress 0.2265481
## Run 76 stress 0.2158796
## Run 77 stress 0.2196687
## Run 78 stress 0.2205634
## Run 79 stress 0.2249209
## Run 80 stress 0.2253682
## Run 81 stress 0.2147638
## ... New best solution
## ... Procrustes: rmse 0.02019995 max resid 0.1706894
## Run 82 stress 0.2209526
## Run 83 stress 0.2188914
## Run 84 stress 0.2185411
## Run 85 stress 0.2187632
## Run 86 stress 0.2200865
## Run 87 stress 0.2233313
## Run 88 stress 0.2273555
## Run 89 stress 0.2231699
## Run 90 stress 0.2289956
## Run 91 stress 0.216548
## Run 92 stress 0.2252693
## Run 93 stress 0.2323094
## Run 94 stress 0.2161477
## Run 95 stress 0.2155953
## Run 96 stress 0.224344
## Run 97 stress 0.236495
## Run 98 stress 0.2206388
## Run 99 stress 0.2270426
## Run 100 stress 0.2157051
## Run 101 stress 0.2172895
## Run 102 stress 0.2287911
## Run 103 stress 0.2268189
## Run 104 stress 0.2183477
## Run 105 stress 0.2145317
## ... New best solution
## ... Procrustes: rmse 0.01178818 max resid 0.09079011
## Run 106 stress 0.2161303
## Run 107 stress 0.2194825
## Run 108 stress 0.2280974
## Run 109 stress 0.2240568
## Run 110 stress 0.2205414
## Run 111 stress 0.2177552
## Run 112 stress 0.2227853
## Run 113 stress 0.2223119
## Run 114 stress 0.224426
## Run 115 stress 0.2248239
## Run 116 stress 0.2216196
## Run 117 stress 0.2186093
## Run 118 stress 0.2263006
## Run 119 stress 0.2236608
## Run 120 stress 0.2228005
## Run 121 stress 0.2251688
## Run 122 stress 0.2194827
## Run 123 stress 0.2160562
## Run 124 stress 0.2207825
## Run 125 stress 0.234424
## Run 126 stress 0.2217608
## Run 127 stress 0.224361
## Run 128 stress 0.215184
## Run 129 stress 0.2167913
## Run 130 stress 0.2155703
## Run 131 stress 0.2215338
## Run 132 stress 0.2161727
## Run 133 stress 0.2175302
## Run 134 stress 0.219476
## Run 135 stress 0.2274296
## Run 136 stress 0.2231364
## Run 137 stress 0.2186249
## Run 138 stress 0.2253854
## Run 139 stress 0.2144394
## ... New best solution
## ... Procrustes: rmse 0.008078074 max resid 0.08629875
## Run 140 stress 0.2235408
## Run 141 stress 0.2332037
## Run 142 stress 0.2309434
## Run 143 stress 0.2220927
## Run 144 stress 0.2247971
## Run 145 stress 0.2151208
## Run 146 stress 0.2221761
## Run 147 stress 0.2234113
## Run 148 stress 0.225497
## Run 149 stress 0.2239924
## Run 150 stress 0.2145903
## ... Procrustes: rmse 0.01038181 max resid 0.09337544
## Run 151 stress 0.224213
## Run 152 stress 0.2266936
## Run 153 stress 0.2222239
## Run 154 stress 0.2246215
## Run 155 stress 0.2297331
## Run 156 stress 0.2223295
## Run 157 stress 0.226072
## Run 158 stress 0.2174392
## Run 159 stress 0.2147968
## ... Procrustes: rmse 0.01026084 max resid 0.09249678
## Run 160 stress 0.2241629
## Run 161 stress 0.2217412
## Run 162 stress 0.2145619
## ... Procrustes: rmse 0.008706978 max resid 0.08607246
## Run 163 stress 0.2256237
## Run 164 stress 0.223733
## Run 165 stress 0.2221946
## Run 166 stress 0.2221341
## Run 167 stress 0.2259662
## Run 168 stress 0.2164148
## Run 169 stress 0.2171011
## Run 170 stress 0.2208268
## Run 171 stress 0.2283241
## Run 172 stress 0.2215728
## Run 173 stress 0.2217803
## Run 174 stress 0.232179
## Run 175 stress 0.2251544
## Run 176 stress 0.2212043
## Run 177 stress 0.2263419
## Run 178 stress 0.21442
## ... New best solution
## ... Procrustes: rmse 0.001703577 max resid 0.0159139
## Run 179 stress 0.2188793
## Run 180 stress 0.2204985
## Run 181 stress 0.2335261
## Run 182 stress 0.2189319
## Run 183 stress 0.2226606
## Run 184 stress 0.2144252
## ... Procrustes: rmse 0.00148575 max resid 0.01558771
## Run 185 stress 0.226184
## Run 186 stress 0.2265209
## Run 187 stress 0.22642
## Run 188 stress 0.2161475
## Run 189 stress 0.22238
## Run 190 stress 0.2270835
## Run 191 stress 0.218626
## Run 192 stress 0.2155601
## Run 193 stress 0.2226199
## Run 194 stress 0.2144192
## ... New best solution
## ... Procrustes: rmse 0.0007192494 max resid 0.008287785
## ... Similar to previous best
## *** Best solution repeated 1 times
##
## Call:
## metaMDS(comm = mat.dis_fish, distance = "bray", k = 2, trymax = 500)
##
## global Multidimensional Scaling using monoMDS
##
## Data: mat.dis_fish
## Distance: bray
##
## Dimensions: 2
## Stress: 0.2144192
## Stress type 1, weak ties
## Best solution was repeated 1 time in 194 tries
## The best solution was from try 194 (random start)
## Scaling: centring, PC rotation, halfchange scaling
## Species: scores missing
## NMDS1 NMDS2 Site
## 1 -0.01430427 -0.10472333 1
## 2 -0.27097428 -0.16902943 1
## 3 0.18940888 -0.33371044 1
## 4 -0.35126583 -0.03394273 1
## 5 0.19887783 -0.12031055 1
## 6 -0.26864930 -0.38585652 1## $vectors
## NMDS1 NMDS2 r2 Pr(>r)
## Acanthurus nigrofuscus -0.89278 0.45050 0.0497 0.033 *
## Caranx melampygus 0.74587 0.66610 0.0025 0.846
## Chaetodon vagabundus 0.98295 0.18387 0.2516 0.001 ***
## Chlorurus spilurus -0.88929 0.45735 0.0058 0.607
## Gomphosus varius -0.35830 -0.93361 0.1819 0.001 ***
## Halichoeres trimaculatus -0.04269 0.99909 0.2409 0.001 ***
## Mulloidichthys flavolineatus 0.16917 0.98559 0.0194 0.205
## Parupeneus multifasciatus -0.99958 0.02883 0.0205 0.224
## Rhinecanthus aculeatus -0.25302 0.96746 0.2328 0.001 ***
## Stegastes sp -0.37836 -0.92566 0.4281 0.001 ***
## Stethojulis bandanensis 0.00125 1.00000 0.2329 0.001 ***
## Thalassoma hardwicke -0.68698 -0.72668 0.0987 0.002 **
## Acanthurus triostegus 0.82689 -0.56236 0.0536 0.038 *
## Chaetodon ephippium -0.64989 -0.76003 0.0372 0.084 .
## Chaetodon lunula -0.39847 -0.91718 0.0905 0.003 **
## Chaetodon ulietensis -0.31336 -0.94963 0.0190 0.253
## Epinephelus merra -0.74968 -0.66179 0.0356 0.074 .
## Zebrasoma scopas -0.67906 -0.73408 0.0462 0.041 *
## Carcharhinus melanopterus 0.61004 -0.79237 0.0368 0.083 .
## Labroides dimidiatus -0.00553 -0.99998 0.0834 0.007 **
## Scarus psittacus 0.07411 0.99725 0.2602 0.001 ***
## Balistapus undulatus -0.70661 -0.70760 0.0958 0.001 ***
## Centropyge flavissima -0.86432 0.50293 0.0128 0.361
## Chaetodon citrinellus -0.95395 0.29996 0.0201 0.231
## Cheilinus chlorourus -0.11859 0.99294 0.0964 0.001 ***
## Cheilinus trilobatus -0.87302 -0.48768 0.0808 0.005 **
## Naso lituratus 0.48500 0.87452 0.0034 0.749
## Cheilio inermis -0.44741 0.89433 0.0549 0.018 *
## Fistularia commersonii -0.12380 -0.99231 0.0176 0.276
## Cephalopholis argus 0.31322 -0.94968 0.0142 0.342
## Chaetodon auriga -0.41540 -0.90964 0.0468 0.031 *
## Chaetodon reticulatus 0.91575 -0.40174 0.0081 0.485
## Cheilinus undulatus 0.91575 -0.40174 0.0081 0.485
## Epibulus insidiator -0.79560 -0.60582 0.0198 0.258
## Siganus spinus -0.07757 0.99699 0.0992 0.003 **
## Chaetodon lunulatus -0.67021 -0.74217 0.0469 0.038 *
## Abudefduf sexfasciatus -0.98281 0.18463 0.0718 0.006 **
## Parupeneus insularis -0.03978 0.99921 0.0052 0.656
## Scarus forsteni 0.00508 -0.99999 0.0045 0.629
## Abudefduf septemfasciatus 0.10859 -0.99409 0.0119 0.374
## Hipposcarus longiceps 0.71779 -0.69626 0.0116 0.347
## Zanclus cornutus -0.70641 0.70780 0.0233 0.157
## Chrysiptera brownriggii -0.79233 -0.61009 0.0152 0.324
## Ostracion meleagris -0.15580 -0.98779 0.0011 0.913
## Pseudocheilinus evanidus -0.93375 -0.35794 0.0012 0.904
## Ellochelon vaigiensis 0.34069 0.94018 0.0016 0.841
## Chaetodon trifascialis -0.80063 0.59915 0.0013 0.904
## Himantura fai 0.99514 0.09850 0.0403 0.070 .
## Scarus oviceps -0.77445 -0.63263 0.0113 0.367
## Acanthurus guttatus 0.93751 -0.34795 0.0127 0.360
## Chaetodon ornatissimus -0.66193 -0.74957 0.0111 0.360
## Oxycheilinus unifasciatus -0.65600 0.75476 0.0130 0.393
## Mulloidichthys vanicolensis -0.16758 0.98586 0.0728 0.016 *
## Aulostomus chinensis -0.26781 0.96347 0.0089 0.499
## Ostracion cubicus -0.14868 0.98888 0.0138 0.304
## Scarus ghobban -0.01733 0.99985 0.0263 0.159
## Ostorhinchus nigrofasciatus 0.22904 0.97342 0.0087 0.518
## Parrotfish Unknonwn -0.72824 0.68532 0.0150 0.299
## Lethrinus olivaceus 0.64738 0.76216 0.0014 0.911
## Parupeneus ciliatus -0.13709 0.99056 0.0193 0.218
## Lutjanus fulvus -0.98294 0.18395 0.0007 0.934
## Chromis viridis -0.77876 0.62733 0.1071 0.002 **
## Coris aygula -0.68498 0.72857 0.0134 0.328
## Stegastes fasciolatus 0.12554 0.99209 0.0242 0.167
## Stegastes albifasciatus -0.99873 -0.05029 0.0042 0.640
## Naso unicornis -0.13541 0.99079 0.0009 0.928
## Parupeneus barberinus 0.24715 0.96898 0.0691 0.019 *
## Dascyllus aruanus -0.33713 -0.94146 0.0737 0.012 *
## Kyphosus sp. 0.09170 -0.99579 0.0413 0.061 .
## Aluterus scriptus -0.08830 -0.99609 0.0062 0.582
## Arothron meleagris 0.19924 -0.97995 0.0205 0.204
## Cantherhines dumerilii 0.31231 -0.94998 0.0270 0.149
## Dascyllus trimaculatus -0.15452 0.98799 0.0002 0.992
## Labroides bicolor -0.28605 -0.95822 0.0097 0.377
## Crenimugil crenilabrus 0.73033 0.68309 0.1067 0.006 **
## Abudefduf sordidus -0.36435 0.93126 0.0121 0.361
## Bodianus perditio -0.13569 0.99075 0.0219 0.154
## Naso annulatus 0.26227 0.96499 0.0096 0.390
## Thalassoma purpureum 0.13288 -0.99113 0.0056 0.573
## Cheilodipterus quinquelineatus 0.24897 0.96851 0.0703 0.019 *
## Pomacentrus pavo 0.28862 0.95744 0.0487 0.049 *
## Canthigaster bennetti 0.25625 0.96661 0.0350 0.098 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Permutation: free
## Number of permutations: 999
##
## $factors
## NULL
##
## $na.action
## function (object, ...)
## UseMethod("na.action")
## <bytecode: 0x143200d00>
## <environment: namespace:stats>## NMDS1 NMDS2 Species pval
## Acanthurus nigrofuscus -0.19901749 0.10042595 Acanthurus nigrofuscus 0.033
## Chaetodon vagabundus 0.49300770 0.09222321 Chaetodon vagabundus 0.001
## Gomphosus varius -0.15281520 -0.39818678 Gomphosus varius 0.001
## Halichoeres trimaculatus -0.02095273 0.49037691 Halichoeres trimaculatus 0.001
## Rhinecanthus aculeatus -0.12207367 0.46676169 Rhinecanthus aculeatus 0.001
## Stegastes sp -0.24756748 -0.60567878 Stegastes sp 0.001
## NMDS1 NMDS2 Species
## Chaetodon vagabundus 0.4930077027 0.09222321 Chaetodon vagabundus
## Gomphosus varius -0.1528151985 -0.39818678 Gomphosus varius
## Halichoeres trimaculatus -0.0209527303 0.49037691 Halichoeres trimaculatus
## Rhinecanthus aculeatus -0.1220736671 0.46676169 Rhinecanthus aculeatus
## Stegastes sp -0.2475674783 -0.60567878 Stegastes sp
## Stethojulis bandanensis 0.0006038459 0.48255616 Stethojulis bandanensis
## Thalassoma hardwicke -0.2158322637 -0.22830545 Thalassoma hardwicke
## Chaetodon lunula -0.1198407201 -0.27584401 Chaetodon lunula
## Labroides dimidiatus -0.0015975728 -0.28870745 Labroides dimidiatus
## Scarus psittacus 0.0378040041 0.50870037 Scarus psittacus
## Balistapus undulatus -0.2186643129 -0.21897099 Balistapus undulatus
## Cheilinus chlorourus -0.0368240043 0.30831429 Cheilinus chlorourus
## Cheilinus trilobatus -0.2481364886 -0.13861335 Cheilinus trilobatus
## Siganus spinus -0.0244356989 0.31406923 Siganus spinus
## Abudefduf sexfasciatus -0.2633377509 0.04947046 Abudefduf sexfasciatus
## Chromis viridis -0.2548157662 0.20526704 Chromis viridis
## Crenimugil crenilabrus 0.2385468294 0.22311580 Crenimugil crenilabrus
## pval
## Chaetodon vagabundus 0.001
## Gomphosus varius 0.001
## Halichoeres trimaculatus 0.001
## Rhinecanthus aculeatus 0.001
## Stegastes sp 0.001
## Stethojulis bandanensis 0.001
## Thalassoma hardwicke 0.002
## Chaetodon lunula 0.003
## Labroides dimidiatus 0.007
## Scarus psittacus 0.001
## Balistapus undulatus 0.001
## Cheilinus chlorourus 0.001
## Cheilinus trilobatus 0.005
## Siganus spinus 0.003
## Abudefduf sexfasciatus 0.006
## Chromis viridis 0.002
## Crenimugil crenilabrus 0.006
Constrained Ordination - Distance based redundancy analysis
I started with using the pulse pc axis, the low tide mean silicate, the distance to seep, and the distance to shore:
Tests revealed the shore distance didn’t matter (permutation test p = 0.555), so I changed this to only use the pulse pc axis, the low tide mean silicate, the distance to seep, and the distance to shore:
## Call: capscale(formula = transformed_communities ~
## Low_Tide_Mean_Phosphate_umolL + Low_Tide_Mean_Silicate_umolL +
## dist_to_seep_m + Low_Tide_Mean_NN_umolL + Low_Tide_Mean_Ammonia_umolL,
## data = transformed_communities_with_explan, distance = "bray")
##
## Inertia Proportion Rank
## Total 32.5200
## RealTotal 41.1784 1.0000
## Constrained 5.9953 0.1456 5
## Unconstrained 35.1832 0.8544 58
## Imaginary -8.6584
## Inertia is squared Bray distance
## Species scores projected from 'transformed_communities'
##
## Eigenvalues for constrained axes:
## CAP1 CAP2 CAP3 CAP4 CAP5
## 2.8690 1.7850 0.6992 0.3968 0.2453
##
## Eigenvalues for unconstrained axes:
## MDS1 MDS2 MDS3 MDS4 MDS5 MDS6 MDS7 MDS8
## 4.145 3.466 2.461 2.069 1.896 1.814 1.432 1.313
## (Showing 8 of 58 unconstrained eigenvalues)The output reports the total inertia, which is the total amount of variation (dissimilarity) in the data. This inertia is decomposed into ‘constrained’ and ‘unconstrained’ components. The constrained component is the total amount of variation explained by the predictors (18.36%), while the unconstrained component is the remaining ‘residual’ variation. There is also info on ‘real’ and ‘imaginary’ components, due to the negative eigenvalues issue which arises with PCoA.
The default plot shows how these variables are loaded onto the first
two CAP axes, and shows how the samples (circles) are ordinated on those
axes, as well as the species scores (red crosses).
We can get the variance explained by these axes from summary:
## $importance
## Importance of components:
## CAP1 CAP2 CAP3 CAP4 CAP5 MDS1 MDS2
## Eigenvalue 2.86899 1.78496 0.69922 0.396820 0.245292 4.1455 3.46634
## Proportion Explained 0.06967 0.04335 0.01698 0.009637 0.005957 0.1007 0.08418
## Cumulative Proportion 0.06967 0.11302 0.13000 0.139636 0.145593 0.2463 0.33044
## MDS3 MDS4 MDS5 MDS6 MDS7 MDS8 MDS9
## Eigenvalue 2.46060 2.06927 1.89624 1.81364 1.43191 1.31320 1.1528
## Proportion Explained 0.05975 0.05025 0.04605 0.04404 0.03477 0.03189 0.0280
## Cumulative Proportion 0.39020 0.44045 0.48650 0.53054 0.56531 0.59721 0.6252
## MDS10 MDS11 MDS12 MDS13 MDS14 MDS15 MDS16
## Eigenvalue 1.08394 0.99752 0.93176 0.88247 0.81741 0.77696 0.7126
## Proportion Explained 0.02632 0.02422 0.02263 0.02143 0.01985 0.01887 0.0173
## Cumulative Proportion 0.65152 0.67575 0.69838 0.71981 0.73966 0.75852 0.7758
## MDS17 MDS18 MDS19 MDS20 MDS21 MDS22 MDS23
## Eigenvalue 0.68219 0.66369 0.61093 0.57537 0.54808 0.51567 0.47836
## Proportion Explained 0.01657 0.01612 0.01484 0.01397 0.01331 0.01252 0.01162
## Cumulative Proportion 0.79240 0.80851 0.82335 0.83732 0.85063 0.86316 0.87477
## MDS24 MDS25 MDS26 MDS27 MDS28 MDS29
## Eigenvalue 0.44304 0.394182 0.356226 0.333264 0.31214 0.293873
## Proportion Explained 0.01076 0.009573 0.008651 0.008093 0.00758 0.007137
## Cumulative Proportion 0.88553 0.895103 0.903754 0.911847 0.91943 0.926564
## MDS30 MDS31 MDS32 MDS33 MDS34 MDS35
## Eigenvalue 0.271585 0.251402 0.243443 0.222833 0.207631 0.190622
## Proportion Explained 0.006595 0.006105 0.005912 0.005411 0.005042 0.004629
## Cumulative Proportion 0.933159 0.939265 0.945176 0.950588 0.955630 0.960259
## MDS36 MDS37 MDS38 MDS39 MDS40 MDS41
## Eigenvalue 0.184730 0.169979 0.149546 0.13959 0.117512 0.110795
## Proportion Explained 0.004486 0.004128 0.003632 0.00339 0.002854 0.002691
## Cumulative Proportion 0.964745 0.968873 0.972505 0.97589 0.978748 0.981439
## MDS42 MDS43 MDS44 MDS45 MDS46 MDS47
## Eigenvalue 0.106832 0.098982 0.082899 0.073070 0.06464 0.059520
## Proportion Explained 0.002594 0.002404 0.002013 0.001774 0.00157 0.001445
## Cumulative Proportion 0.984033 0.986437 0.988450 0.990225 0.99179 0.993240
## MDS48 MDS49 MDS50 MDS51 MDS52 MDS53
## Eigenvalue 0.053156 0.048698 0.0383258 0.035949 0.0265654 0.0205137
## Proportion Explained 0.001291 0.001183 0.0009307 0.000873 0.0006451 0.0004982
## Cumulative Proportion 0.994531 0.995714 0.9966443 0.997517 0.9981624 0.9986606
## MDS54 MDS55 MDS56 MDS57 MDS58
## Eigenvalue 0.0177533 0.0143017 0.0126491 0.0076031 2.847e-03
## Proportion Explained 0.0004311 0.0003473 0.0003072 0.0001846 6.915e-05
## Cumulative Proportion 0.9990917 0.9994390 0.9997462 0.9999308 1.000e+00The first CAP axis explains 6.97% of the constrained community variation, and the second axis explains 4.34%. Therefore, this represents 2.07% of the total community variation.
Though this is already pretty poor, I might as well finish it out as a coding exercise. The loading coefficients for the explanatory variables on the constrained axes are:
## CAP1 CAP2 CAP3 CAP4
## Low_Tide_Mean_Phosphate_umolL 0.1378526 0.10287882 0.0236568 0.85488628
## Low_Tide_Mean_Silicate_umolL 0.1348215 -0.07396710 -0.4770118 0.84007095
## dist_to_seep_m -0.9226700 0.22245994 -0.1614185 -0.04310602
## Low_Tide_Mean_NN_umolL -0.3921606 0.09947864 -0.3572339 0.82416228
## Low_Tide_Mean_Ammonia_umolL -0.2603603 -0.59623153 0.3245798 0.52529160
## CAP5 MDS1
## Low_Tide_Mean_Phosphate_umolL 0.4888992 0
## Low_Tide_Mean_Silicate_umolL 0.2075876 0
## dist_to_seep_m 0.2669787 0
## Low_Tide_Mean_NN_umolL -0.1716233 0
## Low_Tide_Mean_Ammonia_umolL 0.4420828 0To test the important of all the predictors in combination:
## Permutation test for capscale under reduced model
## Permutation: free
## Number of permutations: 999
##
## Model: capscale(formula = transformed_communities ~ Low_Tide_Mean_Phosphate_umolL + Low_Tide_Mean_Silicate_umolL + dist_to_seep_m + Low_Tide_Mean_NN_umolL + Low_Tide_Mean_Ammonia_umolL, data = transformed_communities_with_explan, distance = "bray")
## Df SumOfSqs F Pr(>F)
## Model 5 5.995 4.7031 0.001 ***
## Residual 138 35.183
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1So this tests whether the total variation explained by the constrained axes is significant. Interestingly, though the amount of variation seems trivial, it is showing up as significant…
## Permutation test for capscale under reduced model
## Forward tests for axes
## Permutation: free
## Number of permutations: 999
##
## Model: capscale(formula = transformed_communities ~ Low_Tide_Mean_Phosphate_umolL + Low_Tide_Mean_Silicate_umolL + dist_to_seep_m + Low_Tide_Mean_NN_umolL + Low_Tide_Mean_Ammonia_umolL, data = transformed_communities_with_explan, distance = "bray")
## Df SumOfSqs F Pr(>F)
## CAP1 1 2.869 11.2531 0.001 ***
## CAP2 1 1.785 7.0012 0.001 ***
## CAP3 1 0.699 2.7426 0.002 **
## CAP4 1 0.397 1.5565 0.144
## CAP5 1 0.245 0.9621 0.462
## Residual 138 35.183
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1So, both the axes we have plotted are considered significant.
## Permutation test for capscale under reduced model
## Marginal effects of terms
## Permutation: free
## Number of permutations: 999
##
## Model: capscale(formula = transformed_communities ~ Low_Tide_Mean_Phosphate_umolL + Low_Tide_Mean_Silicate_umolL + dist_to_seep_m + Low_Tide_Mean_NN_umolL + Low_Tide_Mean_Ammonia_umolL, data = transformed_communities_with_explan, distance = "bray")
## Df SumOfSqs F Pr(>F)
## Low_Tide_Mean_Phosphate_umolL 1 1.225 4.8064 0.001 ***
## Low_Tide_Mean_Silicate_umolL 1 0.856 3.3578 0.001 ***
## dist_to_seep_m 1 1.052 4.1245 0.001 ***
## Low_Tide_Mean_NN_umolL 1 0.557 2.1854 0.007 **
## Low_Tide_Mean_Ammonia_umolL 1 1.743 6.8366 0.001 ***
## Residual 138 35.183
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Interestingly, all the axes are coming up as significant contributors to this variation.
PERMANOVA:
## Permutation test for adonis under reduced model
## Marginal effects of terms
## Permutation: free
## Number of permutations: 999
##
## adonis2(formula = transformed_communities ~ Low_Tide_Mean_Phosphate_umolL + Low_Tide_Mean_Silicate_umolL + dist_to_seep_m + Low_Tide_Mean_NN_umolL + Low_Tide_Mean_Ammonia_umolL, data = transformed_communities_with_explan, by = "margin", dist = "bray")
## Df SumOfSqs R2 F Pr(>F)
## Low_Tide_Mean_Phosphate_umolL 1 1.176 0.03616 6.0675 0.001 ***
## Low_Tide_Mean_Silicate_umolL 1 0.807 0.02481 4.1628 0.001 ***
## dist_to_seep_m 1 1.004 0.03088 5.1817 0.001 ***
## Low_Tide_Mean_NN_umolL 1 0.495 0.01522 2.5534 0.006 **
## Low_Tide_Mean_Ammonia_umolL 1 1.701 0.05230 8.7757 0.001 ***
## Residual 138 26.746 0.82244
## Total 143 32.520 1.00000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Again, all of our variables are significant contributors to community composition, though with very low r2.
Mixed Models for Multivariate Data
The first model I test is looking at the probability of presence, and whether it differs between the CV of Salinity. With this, the effect of salinity differs randomly around the species. The random effect for salinity can quantify how much the community composition chnanges along the gradient. If the random effect is non-0, the species change in relative abundance along the salinity gradient.
## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: binomial ( logit )
## Formula: Presence ~ CV_Salinity + (1 + CV_Salinity | species)
## Data: presence_explan
##
## AIC BIC logLik deviance df.resid
## 6771.6 6808.5 -3380.8 6761.6 11803
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.4000 -0.3156 -0.1554 -0.1081 12.7475
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## species (Intercept) 2.859 1.691
## CV_Salinity 1314.560 36.257 1.00
## Number of obs: 11808, groups: species, 82
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.3899 0.2069 -11.55 <2e-16 ***
## CV_Salinity -91.1111 8.8483 -10.30 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## CV_Salinity -0.089## $species
We can test the salinity effect by dropping the random slope and doing an LRT
## Data: presence_explan
## Models:
## mod.no.salinity: Presence ~ CV_Salinity + (1 | species)
## mod.salinity: Presence ~ CV_Salinity + (1 + CV_Salinity | species)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod.no.salinity 3 6769.1 6791.2 -3381.5 6763.1
## mod.salinity 5 6771.6 6808.5 -3380.8 6761.6 1.4393 2 0.4869This doesn’t look great, but the technique is interesting.
GAM
I want to try to do this to look at non-linear distribution along the gradient. I think some species might be non-linear as they might avoid the highly fresh region. I think this might be easiest to look at if I only looked at certain species of interest, which I still would like to define.
For now, I’m setting up the code using species driving differences in the nmds
##
## Family: binomial
## Link function: logit
##
## Formula:
## Presence ~ s(dist_to_seep_m, by = species)
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.06842 0.07563 -27.35 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq
## s(dist_to_seep_m):speciesAbudefduf sexfasciatus 3.232 3.992 4.647
## s(dist_to_seep_m):speciesAcanthurus nigrofuscus 7.603 8.167 116.257
## s(dist_to_seep_m):speciesAcanthurus triostegus 8.464 8.910 39.134
## s(dist_to_seep_m):speciesBalistapus undulatus 8.998 9.000 108.720
## s(dist_to_seep_m):speciesChaetodon auriga 7.494 8.161 59.152
## s(dist_to_seep_m):speciesChaetodon lunula 2.886 3.577 59.245
## s(dist_to_seep_m):speciesChaetodon lunulatus 2.117 2.640 2.369
## s(dist_to_seep_m):speciesChaetodon vagabundus 6.766 7.720 39.815
## s(dist_to_seep_m):speciesCheilinus chlorourus 8.803 8.987 92.577
## s(dist_to_seep_m):speciesCheilinus trilobatus 8.991 9.000 95.455
## s(dist_to_seep_m):speciesCheilio inermis 2.956 3.674 8.693
## s(dist_to_seep_m):speciesCheilodipterus quinquelineatus 1.000 1.001 2.184
## s(dist_to_seep_m):speciesChromis viridis 5.151 6.240 16.583
## s(dist_to_seep_m):speciesCrenimugil crenilabrus 1.760 2.177 9.114
## s(dist_to_seep_m):speciesDascyllus aruanus 1.000 1.001 0.940
## s(dist_to_seep_m):speciesGomphosus varius 6.325 7.217 76.284
## s(dist_to_seep_m):speciesHalichoeres trimaculatus 9.000 9.000 148.689
## s(dist_to_seep_m):speciesLabroides dimidiatus 2.781 3.443 11.877
## s(dist_to_seep_m):speciesMulloidichthys vanicolensis 1.727 2.145 1.660
## s(dist_to_seep_m):speciesParupeneus barberinus 1.000 1.000 2.825
## s(dist_to_seep_m):speciesPomacentrus pavo 1.000 1.000 1.224
## s(dist_to_seep_m):speciesRhinecanthus aculeatus 8.999 9.000 84.167
## s(dist_to_seep_m):speciesScarus psittacus 9.000 9.000 141.910
## s(dist_to_seep_m):speciesSiganus spinus 8.588 8.944 30.934
## s(dist_to_seep_m):speciesStegastes sp 8.387 8.902 15.636
## s(dist_to_seep_m):speciesStethojulis bandanensis 9.000 9.000 157.616
## s(dist_to_seep_m):speciesThalassoma hardwicke 8.647 8.954 76.590
## s(dist_to_seep_m):speciesZebrasoma scopas 6.892 7.898 41.808
## p-value
## s(dist_to_seep_m):speciesAbudefduf sexfasciatus 0.314397
## s(dist_to_seep_m):speciesAcanthurus nigrofuscus < 2e-16 ***
## s(dist_to_seep_m):speciesAcanthurus triostegus 9.81e-06 ***
## s(dist_to_seep_m):speciesBalistapus undulatus < 2e-16 ***
## s(dist_to_seep_m):speciesChaetodon auriga < 2e-16 ***
## s(dist_to_seep_m):speciesChaetodon lunula < 2e-16 ***
## s(dist_to_seep_m):speciesChaetodon lunulatus 0.340470
## s(dist_to_seep_m):speciesChaetodon vagabundus < 2e-16 ***
## s(dist_to_seep_m):speciesCheilinus chlorourus < 2e-16 ***
## s(dist_to_seep_m):speciesCheilinus trilobatus < 2e-16 ***
## s(dist_to_seep_m):speciesCheilio inermis 0.054890 .
## s(dist_to_seep_m):speciesCheilodipterus quinquelineatus 0.139458
## s(dist_to_seep_m):speciesChromis viridis 0.011276 *
## s(dist_to_seep_m):speciesCrenimugil crenilabrus 0.014639 *
## s(dist_to_seep_m):speciesDascyllus aruanus 0.332595
## s(dist_to_seep_m):speciesGomphosus varius < 2e-16 ***
## s(dist_to_seep_m):speciesHalichoeres trimaculatus < 2e-16 ***
## s(dist_to_seep_m):speciesLabroides dimidiatus 0.013746 *
## s(dist_to_seep_m):speciesMulloidichthys vanicolensis 0.473048
## s(dist_to_seep_m):speciesParupeneus barberinus 0.092850 .
## s(dist_to_seep_m):speciesPomacentrus pavo 0.268596
## s(dist_to_seep_m):speciesRhinecanthus aculeatus < 2e-16 ***
## s(dist_to_seep_m):speciesScarus psittacus < 2e-16 ***
## s(dist_to_seep_m):speciesSiganus spinus 0.000299 ***
## s(dist_to_seep_m):speciesStegastes sp 0.071889 .
## s(dist_to_seep_m):speciesStethojulis bandanensis < 2e-16 ***
## s(dist_to_seep_m):speciesThalassoma hardwicke < 2e-16 ***
## s(dist_to_seep_m):speciesZebrasoma scopas < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.38 Deviance explained = 34.2%
## UBRE = -0.11481 Scale est. = 1 n = 4032
##
## Family: binomial
## Link function: logit
##
## Formula:
## Presence ~ s(dist_to_seep_m, by = species)
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -3.058 0.101 -30.29 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq
## s(dist_to_seep_m):speciesAcanthurus guttatus 1.000 1.001 0.342
## s(dist_to_seep_m):speciesAcanthurus nigrofuscus 7.672 8.213 202.150
## s(dist_to_seep_m):speciesAcanthurus triostegus 8.914 8.997 84.671
## s(dist_to_seep_m):speciesCaranx melampygus 8.784 8.985 62.379
## s(dist_to_seep_m):speciesCephalopholis argus 1.000 1.000 1.382
## s(dist_to_seep_m):speciesCheilinus chlorourus 8.919 8.998 169.595
## s(dist_to_seep_m):speciesCheilinus trilobatus 8.999 9.000 189.746
## s(dist_to_seep_m):speciesCheilinus undulatus 1.000 1.000 0.477
## s(dist_to_seep_m):speciesChlorurus spilurus 8.983 9.000 252.182
## s(dist_to_seep_m):speciesEpinephelus merra 5.874 6.793 32.542
## s(dist_to_seep_m):speciesKyphosus sp. 1.000 1.000 1.032
## s(dist_to_seep_m):speciesLethrinus olivaceus 2.545 3.126 12.497
## s(dist_to_seep_m):speciesLutjanus fulvus 1.000 1.000 0.000
## s(dist_to_seep_m):speciesMulloidichthys flavolineatus 8.999 9.000 170.476
## s(dist_to_seep_m):speciesMulloidichthys vanicolensis 2.160 2.664 5.367
## s(dist_to_seep_m):speciesNaso annulatus 1.000 1.000 0.311
## s(dist_to_seep_m):speciesNaso lituratus 3.065 3.760 4.717
## s(dist_to_seep_m):speciesNaso unicornis 1.000 1.000 0.093
## s(dist_to_seep_m):speciesParupeneus barberinus 1.000 1.000 6.975
## s(dist_to_seep_m):speciesParupeneus ciliatus 1.000 1.000 0.025
## s(dist_to_seep_m):speciesParupeneus insularis 1.000 1.000 4.662
## s(dist_to_seep_m):speciesParupeneus multifasciatus 8.364 8.727 157.314
## s(dist_to_seep_m):speciesScarus forsteni 1.000 1.000 0.416
## s(dist_to_seep_m):speciesScarus ghobban 1.000 1.000 0.749
## s(dist_to_seep_m):speciesScarus oviceps 1.000 1.000 0.224
## s(dist_to_seep_m):speciesScarus psittacus 9.000 9.000 250.520
## s(dist_to_seep_m):speciesSiganus spinus 8.941 8.999 93.011
## p-value
## s(dist_to_seep_m):speciesAcanthurus guttatus 0.55925
## s(dist_to_seep_m):speciesAcanthurus nigrofuscus < 2e-16 ***
## s(dist_to_seep_m):speciesAcanthurus triostegus < 2e-16 ***
## s(dist_to_seep_m):speciesCaranx melampygus < 2e-16 ***
## s(dist_to_seep_m):speciesCephalopholis argus 0.23973
## s(dist_to_seep_m):speciesCheilinus chlorourus < 2e-16 ***
## s(dist_to_seep_m):speciesCheilinus trilobatus < 2e-16 ***
## s(dist_to_seep_m):speciesCheilinus undulatus 0.48976
## s(dist_to_seep_m):speciesChlorurus spilurus < 2e-16 ***
## s(dist_to_seep_m):speciesEpinephelus merra 2.99e-05 ***
## s(dist_to_seep_m):speciesKyphosus sp. 0.30976
## s(dist_to_seep_m):speciesLethrinus olivaceus 0.00679 **
## s(dist_to_seep_m):speciesLutjanus fulvus 0.98432
## s(dist_to_seep_m):speciesMulloidichthys flavolineatus < 2e-16 ***
## s(dist_to_seep_m):speciesMulloidichthys vanicolensis 0.10014
## s(dist_to_seep_m):speciesNaso annulatus 0.57715
## s(dist_to_seep_m):speciesNaso lituratus 0.26205
## s(dist_to_seep_m):speciesNaso unicornis 0.76052
## s(dist_to_seep_m):speciesParupeneus barberinus 0.00827 **
## s(dist_to_seep_m):speciesParupeneus ciliatus 0.87343
## s(dist_to_seep_m):speciesParupeneus insularis 0.03085 *
## s(dist_to_seep_m):speciesParupeneus multifasciatus < 2e-16 ***
## s(dist_to_seep_m):speciesScarus forsteni 0.51886
## s(dist_to_seep_m):speciesScarus ghobban 0.38695
## s(dist_to_seep_m):speciesScarus oviceps 0.63580
## s(dist_to_seep_m):speciesScarus psittacus < 2e-16 ***
## s(dist_to_seep_m):speciesSiganus spinus < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.329 Deviance explained = 31.9%
## UBRE = -0.32658 Scale est. = 1 n = 3888
This code didn’t really work as intended. But, when I talked to Kyle Edwards, he suggested that I make a different gam for each species. This is more or less a preliminary step to see if a linear fit is good enough for this data:
From further inspection, the model has issues with Naso annulatus, Naso unicornis, Scarus oviceps
print(warnings_list[[“Naso unicornis”]]) [1] “Iteration limit reached without full convergence - check carefully”
print(warnings_list[[“Naso annulatus”]]) [1] “Fitting terminated with step failure - check results carefully”
print(warnings_list[[“Scarus oviceps”]]) [1] “Fitting terminated with step failure - check results carefully”